ldat <- lapply(ldat,function(x) {
colnames(x) <- gsub(".bam","",gsub(".*:","",colnames(x)))
x
})
WishBone_tsne <- read.csv(file="/Volumes/bioinfomatics$/Mehdi/Matthias/scRNA-seq/WishBone_tsne_with_highTCR.csv", header=TRUE, sep=",")
colnames(WishBone_tsne) <- c('CellName','x','y','Sample')
colData <- read.csv(file="/Volumes/bioinfomatics$/Mehdi/Matthias/scRNA-seq/colData.csv", header=TRUE, sep=",")
merged_meta1 <- merge(x = WishBone_tsne, y = colData, by = "CellName")
color_dt <- data.frame(Sample.x=c('WT1','WT2','WT3','WT4','WT5','highTCR_Rest'), color=c('#f44e42','#ffac38','#bcff37','#36eeff','#9335ff','#36ffee'))
merged_meta <- merge(x = merged_meta1, y = color_dt, by = "Sample.x")
merged_meta <- merged_meta[,c('CellCode','Sample.x')]
rownames(merged_meta) <- merged_meta$CellCode
merged_meta$rep <- 'rep1'
merged_meta <- merged_meta[,c('rep','Sample.x')]
colnames(merged_meta) <- c('rep','celltype')
merged_meta$celltype <- gsub('highTCR_Rest','TCRhiDP',gsub('WT5','CD8SP',gsub("WT4","CD4SP",gsub("WT3","CD4+8low",gsub("WT2","CD69posDP",gsub("WT1","CD69negDP",merged_meta$celltype))))))
# exonic read (spliced) expression matrix
emat_rep1 <- ldat$spliced
colnames(emat_rep1) <- gsub("-","_",colnames(emat_rep1))
emat_rep1 <- emat_rep1[, colnames(emat_rep1) %in% as.vector(rownames(merged_meta))];
# intronic read (unspliced) expression matrix
nmat_rep1 <- ldat$unspliced
colnames(nmat_rep1) <- gsub("-","_",colnames(nmat_rep1))
nmat_rep1 <- nmat_rep1[, colnames(nmat_rep1) %in% as.vector(rownames(merged_meta))];
# spanning read (intron+exon) expression matrix
smat_rep1 <- ldat$spanning;
colnames(smat_rep1) <- gsub("-","_",colnames(smat_rep1))
smat_rep1 <- smat_rep1[, colnames(smat_rep1) %in% as.vector(rownames(merged_meta))];
############
#merging rep1 & rep2
############
emat <- cbind(emat_rep1,emat_rep2)
nmat <- cbind(nmat_rep1,nmat_rep2)
smat <- cbind(smat_rep1,smat_rep2)
meta_info <- rbind(merged_meta,meta_data)
###########
# Running Seurat
###########
WT_Thymocye <- CreateSeuratObject(emat, meta.data = meta_info)
WT_Thymocye.list <- SplitObject(WT_Thymocye, split.by = "rep")
for (i in 1:length(WT_Thymocye.list)) {
WT_Thymocye.list[[i]] <- NormalizeData(WT_Thymocye.list[[i]], verbose = FALSE)
WT_Thymocye.list[[i]] <- FindVariableFeatures(WT_Thymocye.list[[i]], selection.method = "vst",
nfeatures = 2000, verbose = FALSE)
}
Rep1_log <- t(GetAssayData(WT_Thymocye.list[[1]],slot='data',assay = "RNA"))
Rep2_log <- t(GetAssayData(WT_Thymocye.list[[3]],slot='data',assay = "RNA"))
Rep2KO_log <- t(GetAssayData(WT_Thymocye.list[[2]],slot='data',assay = "RNA"))
meta_info$x1 <- 'Others'
meta_info$x2 <- 'Others'
meta_info$x3 <- 'Others'
meta_info$x4 <- 'Others'
meta_info$x5 <- 'Cd8-'
meta_info$x6 <- 'Cd4-'
meta_info$x7 <- 'Itm2a-'
meta_info$x8 <- 'Stat1-'
meta_info$x9 <- 'Others'
meta_info$x10 <- 'Others'
meta_info$x11 <- 'Others'
meta_info$x12 <- 'Others'
colnames(meta_info) <- c("rep", "celltype", "Cd4+Cd8+","Cd4+Cd8-", "Cd4-Cd8+", "Cd4-Cd8-", 'Cd8', 'Cd4', 'Itm2a', 'Stat1', "Zbtb7b+Runx3+","Zbtb7b+Runx3-", "Zbtb7b-Runx3+", "Zbtb7b-Runx3-")
meta_info[rownames(Rep1_log[Rep1_log[,'Zbtb7b'] > 0.1 & Rep1_log[,'Runx3'] > 0.1,]),'Zbtb7b+Runx3+'] <- 'Zbtb7b+Runx3+'
meta_info[rownames(Rep2_log[Rep2_log[,'Zbtb7b'] > 0.1 & Rep2_log[,'Runx3'] > 0.1,]),'Zbtb7b+Runx3+'] <- 'Zbtb7b+Runx3+'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Zbtb7b'] > 0.1 & Rep2KO_log[,'Runx3'] > 0.1,]),'Zbtb7b+Runx3+'] <- 'Zbtb7b+Runx3+'
meta_info[rownames(Rep1_log[Rep1_log[,'Zbtb7b'] > 0.1 & Rep1_log[,'Runx3'] < 0.1,]),'Zbtb7b+Runx3-'] <- 'Zbtb7b+Runx3-'
meta_info[rownames(Rep2_log[Rep2_log[,'Zbtb7b'] > 0.1 & Rep2_log[,'Runx3'] < 0.1,]),'Zbtb7b+Runx3-'] <- 'Zbtb7b+Runx3-'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Zbtb7b'] > 0.1 & Rep2KO_log[,'Runx3'] < 0.1,]),'Zbtb7b+Runx3-'] <- 'Zbtb7b+Runx3-'
meta_info[rownames(Rep1_log[Rep1_log[,'Zbtb7b'] < 0.1 & Rep1_log[,'Runx3'] > 0.1,]),'Zbtb7b-Runx3+'] <- 'Zbtb7b-Runx3+'
meta_info[rownames(Rep2_log[Rep2_log[,'Zbtb7b'] < 0.1 & Rep2_log[,'Runx3'] > 0.1,]),'Zbtb7b-Runx3+'] <- 'Zbtb7b-Runx3+'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Zbtb7b'] < 0.1 & Rep2KO_log[,'Runx3'] > 0.1,]),'Zbtb7b-Runx3+'] <- 'Zbtb7b-Runx3+'
meta_info[rownames(Rep1_log[Rep1_log[,'Zbtb7b'] < 0.1 & Rep1_log[,'Runx3'] < 0.1,]),'Zbtb7b-Runx3-'] <- 'Zbtb7b-Runx3-'
meta_info[rownames(Rep2_log[Rep2_log[,'Zbtb7b'] < 0.1 & Rep2_log[,'Runx3'] < 0.1,]),'Zbtb7b-Runx3-'] <- 'Zbtb7b-Runx3-'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Zbtb7b'] < 0.1 & Rep2KO_log[,'Runx3'] < 0.1,]),'Zbtb7b-Runx3-'] <- 'Zbtb7b-Runx3-'
meta_info[rownames(Rep1_log[Rep1_log[,'Cd4'] > 0.1 & Rep1_log[,'Cd8a'] > 0.1,]),'Cd4+Cd8+'] <- 'Cd4+Cd8+'
meta_info[rownames(Rep2_log[Rep2_log[,'Cd4'] > 0.1 & Rep2_log[,'Cd8a'] > 0.1,]),'Cd4+Cd8+'] <- 'Cd4+Cd8+'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Cd4'] > 0.1 & Rep2KO_log[,'Cd8a'] > 0.1,]),'Cd4+Cd8+'] <- 'Cd4+Cd8+'
meta_info[rownames(Rep1_log[Rep1_log[,'Cd4'] > 0.1 & Rep1_log[,'Cd8a'] < 0.1,]),'Cd4+Cd8-'] <- 'Cd4+Cd8-'
meta_info[rownames(Rep2_log[Rep2_log[,'Cd4'] > 0.1 & Rep2_log[,'Cd8a'] < 0.1,]),'Cd4+Cd8-'] <- 'Cd4+Cd8-'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Cd4'] > 0.1 & Rep2KO_log[,'Cd8a'] < 0.1,]),'Cd4+Cd8-'] <- 'Cd4+Cd8-'
meta_info[rownames(Rep1_log[Rep1_log[,'Cd4'] < 0.1 & Rep1_log[,'Cd8a'] > 0.1,]),'Cd4-Cd8+'] <- 'Cd4-Cd8+'
meta_info[rownames(Rep2_log[Rep2_log[,'Cd4'] < 0.1 & Rep2_log[,'Cd8a'] > 0.1,]),'Cd4-Cd8+'] <- 'Cd4-Cd8+'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Cd4'] < 0.1 & Rep2KO_log[,'Cd8a'] > 0.1,]),'Cd4-Cd8+'] <- 'Cd4-Cd8+'
meta_info[rownames(Rep1_log[Rep1_log[,'Cd4'] < 0.1 & Rep1_log[,'Cd8a'] < 0.1,]),'Cd4-Cd8-'] <- 'Cd4-Cd8-'
meta_info[rownames(Rep2_log[Rep2_log[,'Cd4'] < 0.1 & Rep2_log[,'Cd8a'] < 0.1,]),'Cd4-Cd8-'] <- 'Cd4-Cd8-'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Cd4'] < 0.1 & Rep2KO_log[,'Cd8a'] < 0.1,]),'Cd4-Cd8-'] <- 'Cd4-Cd8-'
meta_info[rownames(Rep1_log[Rep1_log[,'Cd8a'] > 0.1,]),'Cd8'] <- 'Cd8+'
meta_info[rownames(Rep2_log[Rep2_log[,'Cd8a'] > 0.1,]),'Cd8'] <- 'Cd8+'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Cd8a'] > 0.1,]),'Cd8'] <- 'Cd8+'
meta_info[rownames(Rep1_log[Rep1_log[,'Cd4'] > 0.1,]),'Cd4'] <- 'Cd4+'
meta_info[rownames(Rep2_log[Rep2_log[,'Cd4'] > 0.1,]),'Cd4'] <- 'Cd4+'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Cd4'] > 0.1,]),'Cd4'] <- 'Cd4+'
meta_info[rownames(Rep1_log[Rep1_log[,'Itm2a'] > 0.5,]),'Itm2a'] <- 'Itm2a+'
meta_info[rownames(Rep2_log[Rep2_log[,'Itm2a'] > 0.5,]),'Itm2a'] <- 'Itm2a+'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Itm2a'] > 0.1,]),'Itm2a'] <- 'Itm2a+'
meta_info[rownames(Rep1_log[Rep1_log[,'Stat1'] > 0.2,]),'Stat1'] <- 'Stat1+'
meta_info[rownames(Rep2_log[Rep2_log[,'Stat1'] > 0.2,]),'Stat1'] <- 'Stat1+'
meta_info[rownames(Rep2KO_log[Rep2KO_log[,'Stat1'] > 0.1,]),'Stat1'] <- 'Stat1+'
meta_info[meta_info$celltype=='CD4SP' & meta_info$Itm2a=='Itm2a+' & meta_info$Stat1=='Stat1-', 'celltype'] <- 'CD4SP_Immature'
meta_info[meta_info$celltype=='CD4SP' & (meta_info$Itm2a=='Itm2a-' | meta_info$Stat1=='Stat1+'), 'celltype'] <- 'CD4SP_Mature'
WT_Thymocye <- CreateSeuratObject(emat, meta.data = meta_info)
WT_Thymocye.list <- SplitObject(WT_Thymocye, split.by = "rep")
for (i in 1:length(WT_Thymocye.list)) {
WT_Thymocye.list[[i]] <- NormalizeData(WT_Thymocye.list[[i]], verbose = FALSE)
WT_Thymocye.list[[i]] <- FindVariableFeatures(WT_Thymocye.list[[i]], selection.method = "vst",
nfeatures = 2000, verbose = FALSE)
}
meta_info <- WT_Thymocye@meta.data
reference.list <- WT_Thymocye.list[c("rep1", "rep2","rep2_KO")]
WT_Thymocye.anchors <- FindIntegrationAnchors(object.list = reference.list, dims = 1:30)
WT_Thymocye.integrated <- IntegrateData(anchorset = WT_Thymocye.anchors, dims = 1:30)
DefaultAssay(WT_Thymocye.integrated) <- "integrated"
# Run the standard workflow for visualization and clustering
WT_Thymocye.integrated <- ScaleData(WT_Thymocye.integrated, verbose = FALSE)
WT_Thymocye.integrated <- RunPCA(WT_Thymocye.integrated, npcs = 30, verbose = FALSE)
WT_Thymocye.integrated <- RunUMAP(WT_Thymocye.integrated, reduction = "pca", dims = 1:30)
gene_set <- read.csv(file="/Volumes/bioinfomatics$/Mehdi/Matthias/scRNA-seq/scRNAseq_CRG_02012020/gene_set/DE_gene_drodeRes_rep2_fset.csv", header=FALSE, sep=",")
emat <- cbind(emat_rep1,emat_rep2)
nmat <- cbind(nmat_rep1,nmat_rep2)
smat <- cbind(smat_rep1,smat_rep2)
# meta_info <- rbind(merged_meta,meta_data)
emat <- emat[intersect(gene_set$V1,rownames(emat)),]
emat <- emat[!(rownames(emat) =='Cd4' | rownames(emat) == 'Cd8a' | rownames(emat) =='Cd8b1'),]
###########
#Running Seurat
###########
WT_Thymocye <- CreateSeuratObject(emat, meta.data = meta_info)
WT_Thymocye.list <- SplitObject(WT_Thymocye, split.by = "rep")
for (i in 1:length(WT_Thymocye.list)) {
WT_Thymocye.list[[i]] <- NormalizeData(WT_Thymocye.list[[i]], verbose = FALSE)
WT_Thymocye.list[[i]] <- FindVariableFeatures(WT_Thymocye.list[[i]], selection.method = "vst",
nfeatures = 2000, verbose = FALSE)
}
reference.list <- WT_Thymocye.list[c("rep1", "rep2","rep2_KO")]
WT_Thymocye.anchors <- FindIntegrationAnchors(object.list = reference.list, dims = 1:30)
WT_Thymocye.integrated <- IntegrateData(anchorset = WT_Thymocye.anchors, dims = 1:30)
# switch to integrated assay. The variable features of this assay are automatically
# set during IntegrateData
DefaultAssay(WT_Thymocye.integrated) <- "integrated"
# Run the standard workflow for visualization and clustering
WT_Thymocye.integrated <- ScaleData(WT_Thymocye.integrated, verbose = FALSE)
WT_Thymocye.integrated <- RunPCA(WT_Thymocye.integrated, npcs = 30, verbose = FALSE)
WT_Thymocye.integrated <- RunUMAP(WT_Thymocye.integrated, reduction = "pca", dims = 1:30)
# Filtering cells with too many or too few cells
Rep1_filtered_rep2 <- subset(WT_Thymocye.integrated, cells = rownames(meta_info[meta_info$rep=='rep1' |
(meta_info$rep=='rep2_KO' & meta_info$celltype=='CD69posDP' & meta_info$nCount_RNA < 1500000 & meta_info$nCount_RNA > 750000) |
(meta_info$rep=='rep2_KO' & meta_info$celltype!='CD69posDP' & meta_info$nCount_RNA < 1500000 & meta_info$nCount_RNA > 500000) |
(meta_info$rep=='rep2' & meta_info$celltype=='CD69posDP' & meta_info$nCount_RNA < 1500000 & meta_info$nCount_RNA > 750000) |
(meta_info$rep=='rep2' & meta_info$celltype!='CD69posDP' & meta_info$nCount_RNA < 1500000 & meta_info$nCount_RNA > 500000)
,]))
# Dividing CD4SP to CD4SP Immature and Mature
meta_info_new <- Rep1_filtered_rep2@meta.data
PCA <- as.data.frame(Rep1_filtered_rep2@reductions$pca@cell.embeddings)
rownames(PCA) <- rownames(Rep1_filtered_rep2@reductions$pca@cell.embeddings)
meta_info_new$PC_1 <- PCA$PC_1
meta_info_new$PC_2 <- PCA$PC_2
meta_info_new$PC_3 <- PCA$PC_3
meta_info_new[meta_info_new$celltype=='CD4+8low' & meta_info_new$PC_2 < meta_info_new$PC_1,'celltype'] <- 'CD4SP_Immature'
meta_info[rownames(meta_info_new),'celltype'] <- meta_info_new$celltype
gene_set <- read.csv(file="/Volumes/bioinfomatics$/Mehdi/Matthias/scRNA-seq/scRNAseq_CRG_02012020/gene_set/DE_gene_drodeRes_rep2_fset.csv", header=FALSE, sep=",")
emat <- cbind(emat_rep1,emat_rep2)
nmat <- cbind(nmat_rep1,nmat_rep2)
smat <- cbind(smat_rep1,smat_rep2)
emat <- emat[intersect(gene_set$V1,rownames(emat)),]
emat <- emat[!(rownames(emat) =='Cd4' | rownames(emat) == 'Cd8a' | rownames(emat) =='Cd8b1'),]
###########
#Running Seurat
###########
WT_Thymocye <- CreateSeuratObject(emat, meta.data = meta_info)
WT_Thymocye.list <- SplitObject(WT_Thymocye, split.by = "rep")
for (i in 1:length(WT_Thymocye.list)) {
WT_Thymocye.list[[i]] <- NormalizeData(WT_Thymocye.list[[i]], verbose = FALSE)
WT_Thymocye.list[[i]] <- FindVariableFeatures(WT_Thymocye.list[[i]], selection.method = "vst",
nfeatures = 2000, verbose = FALSE)
}
reference.list <- WT_Thymocye.list[c("rep1", "rep2", "rep2_KO")]
WT_Thymocye.anchors <- FindIntegrationAnchors(object.list = reference.list, dims = 1:30)
WT_Thymocye.integrated <- IntegrateData(anchorset = WT_Thymocye.anchors, dims = 1:30)
library(ggplot2)
library(cowplot)
# switch to integrated assay. The variable features of this assay are automatically
# set during IntegrateData
DefaultAssay(WT_Thymocye.integrated) <- "integrated"
# Run the standard workflow for visualization and clustering
WT_Thymocye.integrated <- ScaleData(WT_Thymocye.integrated, verbose = FALSE)
WT_Thymocye.integrated <- RunPCA(WT_Thymocye.integrated, npcs = 30, verbose = FALSE)
WT_Thymocye.integrated <- RunUMAP(WT_Thymocye.integrated, reduction = "pca", dims = 1:30)
PCA <- as.data.frame(WT_Thymocye.integrated@reductions$pca@cell.embeddings)
# Original runing
emat <- cbind(emat_rep1,emat_rep2)
nmat <- cbind(nmat_rep1,nmat_rep2)
smat <- cbind(smat_rep1,smat_rep2)
###########
#Running Seurat
###########
WT_Thymocye <- CreateSeuratObject(emat, meta.data = meta_info)
WT_Thymocye.list <- SplitObject(WT_Thymocye, split.by = "rep")
for (i in 1:length(WT_Thymocye.list)) {
WT_Thymocye.list[[i]] <- NormalizeData(WT_Thymocye.list[[i]], verbose = FALSE)
WT_Thymocye.list[[i]] <- FindVariableFeatures(WT_Thymocye.list[[i]], selection.method = "vst",
nfeatures = 5000, verbose = FALSE)
}
reference.list <- WT_Thymocye.list[c("rep1", "rep2", "rep2_KO")]
WT_Thymocye.anchors <- FindIntegrationAnchors(object.list = reference.list, dims = 1:30)
WT_Thymocye.integrated <- IntegrateData(anchorset = WT_Thymocye.anchors, dims = 1:30)
# switch to integrated assay. The variable features of this assay are automatically
# set during IntegrateData
DefaultAssay(WT_Thymocye.integrated) <- "integrated"
WT_Thymocye.integrated <- ScaleData(WT_Thymocye.integrated, verbose = FALSE)
WT_Thymocye.integrated <- RunPCA(WT_Thymocye.integrated, npcs = 30, verbose = FALSE)
WT_Thymocye.integrated <- RunUMAP(WT_Thymocye.integrated, reduction = "pca", dims = 1:30)
WT_Thymocye.integrated@reductions$pca@cell.embeddings <- as.matrix(PCA)
draw_Plots <- function(data_scatter_new, meta_data_f, gene_name){
meta_data_f <- cbind(meta_data_f, as.data.frame(data_scatter_new[,gene_name]))
meta_data_f <- meta_data_f[c(-8,-9)]
colnames(meta_data_f)[13] <- gene_name
meta_data_f$celltype <- factor(meta_data_f$celltype, levels=c("CD8SP", "CD4SP_Mature", "Sel_Intermed_Cd4-Cd8+", "Sel_Intermed_Cd4+Cd8-", "Sel_Intermed_Cd4+Cd8+", "CD69negDP"))
if(gene_name=="H2-K1"){
colnames(meta_data_f)[13] <- "H2K1"
gene_name <- "H2K1"}
if(gene_name=="H2-D1"){
colnames(meta_data_f)[13] <- "H2D1"
gene_name <- "H2D1"}
P1 <- ggbarplot(meta_data_f, x = "celltype", y = gene_name,
add = "mean_se",
fill = "celltype", palette = c("blue", "green", "red", "darkturquoise", "black", "darkgray"),
position = position_dodge(0.8))+
theme_classic()+ylab(paste0(gene_name," Expression"))+
theme(text = element_text(size = 16, face = "bold"))+
theme(axis.text = element_text(size = 16, face = "bold"))+
theme(plot.title = element_text(size = 16,face = "bold"))+ coord_flip()
print(P1)
}
####
Rep1_filtered_rep2 <- subset(WT_Thymocye.integrated, cells = rownames(meta_info[meta_info$celltype != "CD4SP_Immature" &
((meta_info$rep=='rep2' & meta_info$celltype=='CD69posDP' & meta_info$nCount_RNA < 1500000 & meta_info$nCount_RNA > 750000) |
(meta_info$rep=='rep2' & meta_info$celltype!='CD69posDP' & meta_info$nCount_RNA < 1500000 & meta_info$nCount_RNA > 500000))
,]))
Rep1_filtered_rep2 <- subset(Rep1_filtered_rep2, cells = rownames(meta_info[
!((meta_info$celltype=="CD69posDP" | meta_info$celltype=="TCRhiDP" | meta_info$celltype=="CD4+8low") & meta_info$Cd4=="Cd4-" & meta_info$Cd8=="Cd8-")
,]))
Rep1_filtered_rep2@meta.data[((Rep1_filtered_rep2@meta.data$celltype=="CD69posDP" | Rep1_filtered_rep2@meta.data$celltype=="TCRhiDP" | Rep1_filtered_rep2@meta.data$celltype=="CD4+8low") & Rep1_filtered_rep2@meta.data$Cd4=="Cd4+" & Rep1_filtered_rep2@meta.data$Cd8=="Cd8+"), 'celltype'] <- "Sel_Intermed_Cd4+Cd8+"
Rep1_filtered_rep2@meta.data[((Rep1_filtered_rep2@meta.data$celltype=="CD69posDP" | Rep1_filtered_rep2@meta.data$celltype=="TCRhiDP" | Rep1_filtered_rep2@meta.data$celltype=="CD4+8low") & Rep1_filtered_rep2@meta.data$Cd4=="Cd4+" & Rep1_filtered_rep2@meta.data$Cd8=="Cd8-"), 'celltype'] <- "Sel_Intermed_Cd4+Cd8-"
Rep1_filtered_rep2@meta.data[((Rep1_filtered_rep2@meta.data$celltype=="CD69posDP" | Rep1_filtered_rep2@meta.data$celltype=="TCRhiDP" | Rep1_filtered_rep2@meta.data$celltype=="CD4+8low") & Rep1_filtered_rep2@meta.data$Cd4=="Cd4-" & Rep1_filtered_rep2@meta.data$Cd8=="Cd8+"), 'celltype'] <- "Sel_Intermed_Cd4-Cd8+"
DefaultAssay(Rep1_filtered_rep2) <- "RNA"
data_scatter <- as.data.frame(t(GetAssayData(Rep1_filtered_rep2,slot='data',assay = "RNA")))
meta_data_new <- Rep1_filtered_rep2@meta.data
meta_data_new$PC_1 <- 0
meta_data_new$PC_1 <- Rep1_filtered_rep2@reductions$pca@cell.embeddings[rownames(meta_data_new),'PC_1']
meta_data_new$Cd4Cd8 <- "NA"
meta_data_new$Cd4Cd8 <- paste0(meta_data_new$Cd4, meta_data_new$Cd8)
meta_data_new$PC1_range <- "-10_to_-5"
meta_data_new[meta_data_new$PC_1 > -5 & meta_data_new$PC_1 < -2.5, "PC1_range" ] <- "-5_to_-2.5"
meta_data_new[meta_data_new$PC_1 > -2.5 & meta_data_new$PC_1 < 0, "PC1_range" ] <- "-2.5_to_0"
meta_data_new[meta_data_new$PC_1 > 0 & meta_data_new$PC_1 < 2.5, "PC1_range" ] <- "0_to_2.5"
meta_data_new[meta_data_new$PC_1 > 2.5 & meta_data_new$PC_1 < 7.5, "PC1_range" ] <- "2.5_to_7.5"